Evanescent wave approach to diffractive phenomena in convex bilhards with corners 
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What we are going to call in this paper "diffractive phenomena" in billiards is far from being deeply 
understood. These are sorts of singularities that, for example, some kind of corners introduce in the 
energy eigenfunctions. In this paper we use the well-known scaling quantization procedure to study 
them. We show how the scaling method can be applied to convex billiards with corners, taking into 
account the strong diffraction at them and the techniques needed to solve their Helmholtz equation. 
As an example we study a classically pseudointegrable billiard, the truncated triangle. Then we focus 
our attention on the spectral behavior. A numerical study of the statistical properties of high-lying 
energy levels is carried out. It is found that all computed statistical quantities are roughly described 
by the so-called semi-Poisson statistics, but it is not clear whether the semi-Poisson statistics is the 
correct one in the semiclassical limit. 

PACS numbers: 03.65.Ge, 05.45.Mt 



I. INTRODUCTION 

In this work we calculate very high lying eigenval- 
ues of a billiard system using the so-called scaling 
method 0,0,0]. This method has two main advantages, 
it is formulated on the boundary of the billiard (allowing 
matrices of order fc, the wavenumber) and it avoids ze- 
roes searching algorithms. We have solved the Helmholtz 
equation with Dirichlet boundary conditions on the bil- 
liard boundary. A point worth to mention is that any 
eigenfunction is C°° at domain points. At straight seg- 
ments of the boundary, eigenfunctions are reflected as 
odd functions, so that the result satisfies the Dirichlet 
condition. Though being also C°° at the boundary by 
straight pieces, they are not analytical at a vertex when 
two segments of the boundary join (with inner angle 7r/r, 
where r is a noninteger number). This feature leads to 
what we address as diffractive phenomena. 

The main idea behind the scaling method is that 
the trial functions of the variational problem can be 
parametrised by energy (alternatively, by the wavenum- 
ber k). With this idea in mind these functions are ex- 
panded in a suitable scaling basis (plane waves with dif- 
ferent propagation directions are one example and these 
are the elements we use here). Asking the function to be 
zero at the boundary of the billiard is the same as asking 
its norm there to be zero as well. Then, eigenvalues (and 
eigenfunctions) can be obtained by solving a generalised 
eigenvalue problem that involves the quadratic form as- 
sociated with the norm of the function on the boundary. 

In applying the method to this system, evanescent 
waves were needed. This is because real plane wave so- 
lutions to the Helmholtz equation cannot represent all 
the features of diffraction. It has been shown Q that 
an evanescent plane wave, which oscillates along propa- 
gation direction faster than the wavenumber k, can be 



constructed by means of real plane waves. Nevertheless, 
the corresponding superposition is a singular one, sug- 
gesting the direct use of the evanescent functions in the 
basis. This is the way we have dealt with diffraction. 

For simplicity we restrict our considerations to an angle 
of 37r/4. A particularly suitable example to study the 
effect of this kind of corner on the quantum properties 
is the polygonal billiard shown in Fig. ^ This billiard 
is called the truncated triangle. It has been studied in a 
variety of contexts H E S H H ES IH ■ 
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FIG. 1: The truncated triangle. 

The truncated triangle belongs to the class of ratio- 
nal polygons. That are polygons where all angles aj = 
rrijTT/nj between sides are rationally related to tt, where 
rrij, Uj > are relatively prime integers. The free motion 
inside a rational polygon is integrable if nij = 1 for all 
j, which is the case for rectangles, the equilateral trian- 
gles, the 7r/2, 7r/4, 7r/4-triangles and the 7r/2, 7r/3, 7r/6- 
triangles. All other rational polygons are pseudointe- 
grable 5]. Like in integrable systems, the phase s pac e 
is foliated by two-dimensional invariant surfaces |l2lll3j . 
However, the genus of the surfaces is greater than one 
due to critical corners with > 1 . In the case of the 
truncated triangle, the genus is 17 = 2. Roughly speak- 
ing, the invariant surface is a torus with an additional 
handle. 

Quantum signatures of pseudointegrability can be 
found in the energy eigenfunctions JJ.J and in the sta- 
tistical properties of energy levels The energy lev- 
■dafeEdeatnp^tiiie 'ointegrable systems are correlated in con- 



2 



trast to those of integrable systems which are generically 
weh described by the Poissonian random processes [l5|. 
For example, the nearest-neighbor spacing distribution 
of pseudointegrable systems generically displays level re- 
pulsion [5, resembling the Gaussian orthogonal ensem- 
ble (GOE) of random-matrix theory which describes 
fully chaotic systems with time- reversal symmetry |l7j . 
Significant deviations from GOE have been observed first 
in Refs. 0, • It has been suggested that the spectral 
statistics of pseudointegrable systems is another exam- 
ple of critical or intermediate statistics [l^ l20l | . Criti- 
cal statistics appear in many condensed matter problems 
such as in mesoscopic disordered systems at the critical 
point of the metal-insulator transition |2lj |. in systems 
with a few interacting electrons |22| |. and in incommen- 
surate multiwalled carbon nanotubes p^ . 

In Refs. [l^l20l| it has been proposed to use the semi- 
Poisson (SP) statistics as a reference point for critical 
statistics. The SP statistics is defined by a simple con- 
struction: remove every other level from an ordered Pois- 
son sequence |23, 13 ■ The SP statistics is useful because 
it provides explicit formulas for a number of statistical 
quantities which can be compared to the statistical prop- 
erties of a g iven system. For several pseudointegrable 
systems [ll IH EE El IH IM El it has been con- 
firmed that the SP statistics indeed describes the short- 
range level correlations rather well. Comparing the long- 
range level correlations is numerically a difficult task be- 
cause the statistical properties of rational polygons con- 
verge extremely slowly as energy is increased 28, 29]. 
Fortunately, semiclassical periodic-orbit theory allows to 
compute analytically the long-range level correlations, in 
terms of the level compressibility y, of a few special sys- 
tems, like certain right triangles 2»| and the barrier bil- 
liard ^2^ . In the former case x diifers in general from the 
SP result, whereas in the latter case x is in agreement 
with the SP statistics. Nothing is known analytically 
about the generic case. 

Our numerical analysis will show that the truncated 
triangle is well described by the SP statistics but de- 
viations are not negligible. Our analysis extends that 
reported in Ref. in many respects: (i) The scaling 
method allows us to compute more levels, giving a better 
statistics, (ii) Moreover, high-lying energy levels can be 
computed. This puts us in a position to study the rele- 
vant asymptotic regime, (iii) More statistical quantities 
are computed, (iv) The numerical results are compared 
to the SP statistics. 

The paper is organized as follows: in Section m we 
describe the features that are considered when applying 
the scaling method to the system under investigation. In 
Section UTTl we present the statistical studies carried out 
with the data that has been obtained. Finally, Section 
IIVI is devoted to conclusions. 



II. THE SCALING METHOD APPLIED TO THE 
TRUNCATED TRIANGLE 

We are not going to explain the scaling method here, 
and we address the reader to the Appendix and the given 
references for details. 

An eigenfunction of a billiard can be constructed as 
a plane wave superposition. This can include evanes- 
cent waves, i.e., plane waves with complex wave vectors. 
These types of waves should be present in quantum bil- 
liards and can be associated with diffractive phenomena. 
Discontinuities at the boundaries seem to be strongly re- 
lated to the way evanescent functions must be considered 
to solve the problem. 

Several authors have focused their attention on this 
issue H H S m Eil. For polygonal billiards, in the 
generic case, theory suggests that there is no real plane 
wave superposition that can be an eigenfunction 30] . On 
the other hand, good numerical results using only plane 
waves could be found , but working within the region 
of low energies. Evanescent waves were studied in detail 
by Berry |^ in the context of quantum billiards. In this 
approach, the main reason to do it comes from the idea 
of constructing them by continuation of an external scat- 
tering superposition, containing only real plane waves. 
In principle, this does not seem to be possible. Anyway, 
he showed how evanescent waves can be expressed as the 
singular limit of an angular superposition of real plane 
waves. 

In our system, and working at high energies, we need 
to consider evanescent waves explicitely, because a real 
plane wave representation is singular. In the semiclas- 
sical limit, the only way to obtain eigenfunctions that 
include evanescent waves is by considering them in the 
variational problem In the remaining part of this 
section we are going to show some examples of the waves 
considered and the idea behind their selection. 

Plane wave solutions to the Hclmholtz equation can be 
written as 

'(/'(r) — eyip [ik cos {0 + ia)x + iksin [9 + ia)y\ 

with 6 and a real, using coordinates r = (x, y) in the 
plane. We can express this wave in a slightly different 
fashion by distinguishing the real and imaginary parts in 
the exponent, 

Tpir) = exp [ifc cosh (q;)5:] exp [— fcsinh (a)?;]. 

Here, the propagation direction x = xcosO -f ysmO im- 
plies an angle 9 over the x axis. In this direction the func- 
tion has a wavelength that is given by 27r/(fccoshQ:) < 
27: /k. However, in the orthogonal direction y the func- 
tion is an exponential with coefficient — /csinha. 

An example of these functions applied to our system 
can be found in Fig. |21 where we show an evanescent 
wave corresponding to the truncated triangle of param- 
eter value a = (a/S — l)/2 (Note that since the absolute 
size of the system is irrelevant, we measure length scales 
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as dimensionless quantities) . In this case we chose fc = 20 
in order to have a visible wave. For the higher-lying wave 
numbers {k ~ 1000 or greater) these functions are diffi- 
cult to see in the domain since they decay very fast. Nev- 
ertheless, the slower oscillating functions we have used in 
the latter energy region can be seen over the boundary. 

y 



q=0 q=a 




FIG. 2: Evanescent wave for fc = 20 on the truncated triangle 
domain with parameter value a = (^/5 — l)/2. The propagation 
direction x (given by sin 9 = —0.3 in this case) and the correspond- 
ing decaying direction y are both shown by the arrows. Also, we 
point out values and a of the arclength coordinate q. 

Now that the general ideas related to evanescent waves 
have been exposed, we are going to explain the way we 
selected them, and we are also going to show some ex- 
amples for k — 1000. The first thing to point out is that 
diffraction is generated at the critical corner with angle 
Stt/A that is shown in Fig. ^ So, our sets or famihes of 
evanescent waves are centered at this corner. We consider 
three different kinds of functions, all sharing the previ- 
ous property, but decaying to one or the other "side" 
of the corner (in terms of the arclength boundary coor- 
dinate q) in the first place and also an additional one 
that decays to both sides. We use a family of 12 waves 
whose propagation direction 9 is slightly smaller than 
— 7r/4, i.e., around the right hand side of the boundary 
as seen from Fig. ^ another set of 8 waves with almost 
horizontal propagation and finally one that goes along 
an intermediate direction (for this one we took a value 
of sin 61 = —0.3, as can be seen in Fig. In fact, in 
order to appreciate their contribution, it is convenient to 
look at them on the boundary of the billiard. For this 
reason we show in Fig. |3| three examples of evanescent 
waves that are good representatives of these three fami- 
lies. As already mentioned, the exponential decay is from 
the critical corner, which corresponds to q — a in terms 
of the boundary coordinate. This is the effect of consid- 
ering an angle slightly different from —tt/A angle for the 
propagation direction. The same happens with the sec- 
ond family whose example is directed almost along the 
horizontal (left) segment of the boundary. Finally, we 
show the evanescent wave that decays to both sides of 
the corner. 




0.6 a 0.7 q 0.8 



(b) 




0.3 0.5 a 



(c) 



0.61 a 0.63 



FIG. 3; Three examples of evanescent waves over the billiard 
boundary for k = 1000. Vertical axis corresponds to the real part 
of the evanescent wave in arbitrary units. Horizontal axis corre- 
sponds to the boundary arclength coordinate q (dimensionless). In 
(a) we show a wave decaying to the right (propagation direction 
given by sine = -0.682), in (b) one to the left (sin 6 = -0.025), 
and in (c) to both sides (sin0 = —0.3) from q = a, the position of 
the critical corner. 



This approach to the problem has proven to be very 
efficient. As a matter of fact, studying carefully several 
billiard eigenfunctions that show the greatest norm (er- 
ror) over the boundary (without considering evanescent 
waves in their calculation) , we could check that these are 
the main components of the diffracted field. 

Then we apply the symmetries of our system in order 
to get the right contributions. This is easy to implement 
by using the symmetry operations of the C2v group. 

We have taken only up to 21 evanescent waves in order 
to get our results. This is a small number compared to 
the roughly 1000 real plane waves that are used in the 
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highest-lying energy window we have obtained But 
they are key to lower the error of the eigenvalues. This 
shows that, even though the relevance of this contribu- 
tion goes to zero in the semiclassical limit we need to 
consider it in order to resolve individual states. 



III. SPECTRAL STATISTICS 

We here examine the spectral statistics of the quan- 
tized truncated triangle. We consider the generic case 
where a is an irrational number. In the nongeneric case 
of a being rational, the energy spectrum contains a subset 
of Poisson distributed levels fsl . We focus on the param- 
eter value a = (\/5 — l)/2, the reciprocal of the golden 
mean. For other parameter values of a [{s/b — l)/2-\- 0.2, 
2/7r, and 2/tt + 0.1] similar results have been obtained. 
We use energy windows of length 20 000 in five different 
regimes starting with level number 853, 32 124, 89 607, 
149 879, and 190 356 corresponding to the wave number 
k = 100, 600, 1000, 1300, and 1500. 

To study the local fluctuations in the level sequence 
iJi < £'2 < i?3 < . . . it is necessary to remove the sys- 
tematic global energy dependence of the average density. 
To do so, we "unfold" the spectra in the usual way by 
setting En = N{En)] see, e.g., Ref. [s^. N{E) is the 
smooth part of the integrated density of states, i.e., the 
number of levels up to energy E. We approximate N{E) 
by the generalized Weyl's law including perimeter and 
corner corrections 134 



N{E) 



A 

An' 



-E 



L 

47r 



E + C , 



(1) 



where A ~ a + 1/2 is the area of the billiard, L = 2a + 
2 + \/2 is the perimeter, and C = 11/36 is the corner 
correction. The unfolded spectra {E„} have unit mean 
level spacing and are dimensionless. Henceforth, the tilde 
will be suppressed. 



A. Nearest-neighbor spacing distributions 

The most popular statistical quantity in the field of 
quantum chaos is the nearest-neighbor spacing distribu- 
tion. It is defined as the probability density of the spacing 
s (in units of the mean level spacing) between adjacent 
levels, 



1 " 

P(s) = lim - y S{s - E.,+i + E,) 



(2) 



Clearly, the nearest-neighbor spacing distribution is a 
measure of short-range level correlations. We will com- 
pute its integral, the cumulative spacing distribution 



/(s)= \ P{s')ds' 
Jo 



(3) 



For the Poisson statistics Pp{s) = exp (— s) and /p(s) = 
1 — exp(— s), the GOE is well described by the Wigner 
surmise Pw{s) = (7r/2)s exp (— 7rs^/4) and I^{s) = 1 
exp (— 7rs^/4), and for the SP statistics 



Psp(s) = 4se-^^ Isp{s) - 1 - (2,s + 1) 



-2s 



(4) 



At small s, Pgp(s) exhibits a linear increase from zero 
(level repulsion) similar to the Wigner surmise. At large 
s, PsP (s) has an exponential fall-off as the Poisson statis- 
tics. 

In Fig.^we plot the difference between the cumulative 
spacing distribution to the SP result for the lowest-energy 
window. Good agreement with the SP statistics can be 
observed. Figure[Slshows a magnification, containing also 
the other energy windows. We see increasing deviations 
from SP for k = 600 and k — 1000, and then decreasing 
deviations for k — 1300 and k — 1500. No clear trend to 
SP is visible as one goes to higher energies. 
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FIG. 4: Difference between the cumulative spacing distribution of 
the first energy window and the SP result. 
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FIG. 5: The cumulative spacing distribution of all energy windows. 

Similar small deviations from the SP statistics have 
been observed at the metal-insulator transition in the 
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three-dimensional Anderson model. In that case, an en- 
semble average over certain boundary conditions removes 
the deviations considerably To see whether this is 
possible also in our case, we consider Neumann bound- 
ary conditions on the two boundary segments which do 
not touch the critical corner; see Fig. ^ In this way 
we obtain four energy spectra corresponding to Dirich- 
let/Dirichlet, Dirichlet/Neumann, Neumann/Dirichlet, 
and Neumann/Neumann-boundary conditions. In con- 
trast to the case of the Anderson model we find that av- 
eraging over these boundary conditions does not reduce 
the deviation from the SP statistics. 



B. Next-to-nearest spacing distributions 

We now consider the next-to-nearest spacing distribu- 
tion (second- neighbor-spacin g di stribution) and its inte- 
gral. The SP statistics gives 



Psp(2,s) 
/SP(2,S) 



1 



2s 



{As^ + 6s^ -I- 6s -I- 3)e" 



'2s 



(5) 



Figure El shows that the cumulative next-to-nearest 
spacing distribution is well described by the SP statis- 
tics. The maximal deviations from SP are smaller than 
in the case of the nearest-neighbor spacing distribution 
in Fig.|Sl But, again, there is no clear convergence to the 
SP statistics as one goes to higher energies. 




S(L) is the local variance of the number n{L,E) = 
N{E + L/2) - N{E - L/2) of energy levels in the in- 
terval - L12,E + L/2]. For the SP statistics we 

have [Talis El 



(7) 



Figure shows the number variance computed for the 
five energy windows. Note that the regime is well be- 
low the crossover region where S(L) begins to saturate 
at a value determined by the shortest periodic orbit |36j| . 
With increasing energy the number variance comes closer 
to the SP result, without showing a clear stabilization. 
To estimate the limit curve as fc ^ oo we use the ex- 
trapolation procedure described in Ref. : extrapolate 
point by point (with L fixed) the four highest curves with 
a fit A{L) -\- B{L)/k. The limit curve, A{L), is shown as 
dashed curve in Fig.[7| It is closer to the SP result. To 
measure the difference we concentrate on the slope in the 
region of large L, the so-called level compressibility 



X = lim 

L — >oc 



S(i) 



(8) 



We get X ~ 0-41 from the extrapolated curve. This is 
close to the SP result of 1/2 (1 for Poisson j36|J and 
for GOE [11 HI). 
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semi-Poisson 
k=100-k=1500 
extrapolated 




FIG. 6: The cumulative next-to-nearest spacing distribution. 



FIG. 7: Number variance S(L) for k = 100, 600, 1000, 1300, and 
1500 (from below). 



C. Number variance 

So far we have studied short-range level correla- 
tions (nearest and next-to-nearest spacing distributions) . 
Long-range level correlations are conveniently studied 
with the help of the number variance 

nL) = {{n{L,E)-Lf) . (6) 



D. The form factor 

Another measure of long-range level correlations is the 
form factor K(t), the Fourier transform of the two-point 
correlation function. The limit t ^ is related to 
the number variance by means of K(0) — x (see, e.g., 
Ref. [H) with X from Eq. ©. 
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The form factor can be approximated numerically by 
(see, e.g., Ref. 113) 



K{T;n) 



l+n 
3=1 



(9) 



In our case n = 20 000. Note that r is dimensionless. 
Figure |S1 shows K(T;n) averaged over small intervals of 
size At = 0.006 in the high-energy regime, i.e., k = 1500. 
It is difficult to estimate K(0) from such kinds of noisy 
data, but it is justified to say that K{0) is below the SP 
prediction 1/2, in agreement with our former numerical 
results on the number variance. 

A better way to compare the form factor to the SP 
statistics is introduced in Ref. ,28J. Fit K{T]n) to the 
function 



i^fit(T) = 



- 2c + Att^t^ 



(10) 



If c = 4 then function (|10|) is the SP form factor. We use 
the quantity Xfit(O) — 1/2 to measure the difference to 
the SP statistics. Keep in mind that iffit(O), in general, 
differs from if (0; n) since it depends also on K(t\ n) with 
T > 0. Figure |H1 shows the result obtained by fitting 
Eq. 1)1 0|l to the smoothed data over the range < r < 3. 
We get Ktit{0) ~ 0.552. For the lower-energy windows 
with k = 600, k = 1000, and k = 1300 we find 0.548, 
0.576, and 0.56, respectively. 

Again, as for the other spectral quantities described in 
the previous subsections, we find that the spectral statis- 
tics is roughly described by the SP statistics. However, 
there are small but significant deviations which show no 
clear trend to zero as the energies are increased. 




FIG. 8: The form factor jSJ for the highest -energy window aver- 
aged over small intervals of size At = 0.006. The smooth curve is 
the fit OSJ with c = 4.464. 



IV. CONCLUSION 

We have shown how the scaling method can be applied 
to convex billiards with corners. To gain insight into the 
strong diffractive phenomena that appear in these sys- 
tems we have studied the truncated triangle. As a result 
of our investigations of the diffracted field we could iden- 
tify its main components. These were directly introduced 
in the function basis in order to obtain the spectral data 
of this system. Evanescent waves conveniently selected 
and associated with the centers of diffraction (in our case 
the two straight segments junction at 37r/4) have been 
successfully identified as a very efficient way to deal with 
these kinds of phenomena. This allowed us to obtain a 
great number of highly excited eigenvalues. 

We have studied the statistical properties of high-lying 
energy levels in this pseudointcgrable billiard. We have 
found that the nearest-neighbor spacing distributions, 
next-to-nearest spacing distributions, number variance, 
and the spectral form factor are roughly described by the 
semi-Poisson (SP) statistics. Whether the SP statistics 
is asymptotically the exact statistics cannot be decided. 
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APPENDIX: THE SCALING METHOD 

In this appendix we are going to briefiy explain the 
scaling method Q, |^ |^ . The main point that one should 
have in mind is that the boundary norm can be written as 
a function of energy because scaling is used. If (p{r) sat- 
isfies the Helmholtz equation with eigenvalue fcp and we 
associate with it the scaling function (j>{k,r) = 0(fcr/fco) 
then, these functions satisfy the same equation with 
eigenvalue k^. If we have a billiard defined by a star- 
shaped domain and a A;^ exists such that (ji^ik^^r) — 
at the boundary C, then this is a scaling eigenfunction or, 
equivalently, an eigenfunction of the billiard with Dirich- 
let boundary conditions. 

The boundary norm defined by the expression f{k) = 
§^<l)'^{k,r)dl/rn can be expanded up to third order 
around k^, independently on the exact shape of 0. We re- 
mind that r„ = r • n (where n is the unit outgoing normal 
vector to C) is always positive for star-shaped domains. 
Taking into account this result we can evaluate the norm 
and its derivative in k at value kg = i 



Sfj,, obtaining 



/(fco) 



Sudf 



(11) 
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This useful expansion, dropping terms of order 5^, 
turns out to be our quantization condition. Then, all the 
scaling eigenfunctions with eigenvalues close to ko can be 
found by solving a generalized eigenvalue problem: 



^{ko)-X,F{ko) 



(12) 



In this equation and for numerical calculations, the 
quadratic form F associated with / can be evaluated in a 



basis of scaling functions ■tjji{k, r);i = 1, . . . , N (like plane 
waves for instance) by means of 



Fij{ko) = j) ^i{ko,r)'ipj{ko,r)dl/rn ■ (13) 



The eigenfunctions are (/)^(fc,r) = X^ili V'i(fc7 r), and 
the eigenvalues can be found 
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